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NUMERICAL SPATIAL MARCHING TECHNIQUES FOR ESTIMATING DUCT 
ATTENUATION AND SOURCE PRESSURE PROFILES 
by K. J. Baumeister 
NASA Lc»wis Research Center 
Cleveland, Ohio 44135 
ABSTRACT 

A numerical method is developed that could predict the pressure distribution 
of a ducted source from far- field pressure inputs. Using an initial value formu- 
lation, the two-dimensional homogeneous Helmholtz wave equation (no steady 
flow) is solved using explicit marching techniques. The Von Neumann method 
is used to develop relationships which describe how sound frequency and grid 
spacing effect numerical stability. At the present time, stability considerations 
limit the approach to high frequency sound. Sample calculations for both hard 
and soft wall ducts compare favorably to known boundary value solutions. In 
addition, assuming that reflections in the duct arc small, this initial value ap- 
proach is successfully used to determine the attenuition of a straight soft wall 
duct. Compared to conventional finite difference cr finite element boundary 
value approaches, the numerical marching techni(]ue is orders of magnitude 
shorter in computation time and required computer storage and can be easily 
employed in problems involving high frequency sound. 

INTRODUCTION 

A detailed knowledge of a fan noise source is required in the d;;sign of 
turbojet engine acoustic suppressors, since the energy attenuation of the sup- 
pressor is strongly dependent on the radial and circumferential pressure dis- 
tribution of the noise source (ref. 1). Recent experimental and analytical 
research (refs. 2 to 5) has focused on determining the radial and circumfer- 
ential variations in acoustic pressure for typical turbofan sources from both 
internal and far-field pressure measurements. 

The direct calculation of the internal structure of a ducted noise source 
from far-field pressure measurements is essentially an initial value problem. 
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If pressure and impedance are known at the boundary of the system (far-field), 

« 

the pressure can be uniquely determined in the vicinity of the inlet as well as 
inside the inlet ducting which may contain acoustic energy absorbers in the 
wall. As an alternate to the statistical approach of Ref. 5, a relatively simple 
numerical marching method is developed by which this information, if available 
in the far- field, could be translated directly into a description of a ducted noise 
source. 

In the absence of mean flow, explicit initial value numerical marching 
techniques are developed to obtain a solution to the homogeneous Helmholtz 
equation for a two-dimensional cartesian coordinate system. The technique uses 
the finite difference representation of the Helmholtz equation. Explicit marching 
schemes arc emphasized since they reduce the solution of the Helmholtz equation 
to a simple algebraic algorithm which requires a minimum of computer storage 
and running time. Several sample calculations are also presented to illustrate 
how a source could be predicted from duct far-field conditions. 

In addition, assuming that reflections in the duct are negligible, the technique 
is used to determine the attenuation of a straight two-dimensional soft wall duct 
by prescribing a pressure distribution and the acoustic impedar^ce of the medium 
(p*c*) as the initial conditions at the source (see Symbols for all symbol meanings). 

The truncation error and stability o.*' thetic various marching schemes are 
studied. In ' articular, the Von Neumann method is used to develop relationships 
which describe how sound frequency, transverse grid spacing and marching step 
size affect numerical stability. 

SYMBOLS 

Fourier coefficient. I.q. Bl 
a quadratic constimt, Eq. 26 

a, ^ difference coefficient 

m 

Fourier coefficient, Eq. 18 

b, _ difference coefficient 

m 

c quadratic constant. Eq. 27 
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difference coefficient 
speed of sound 
f|arameter, Eq. 28 
difference coefficient 
sound attenuation 
acoustic power 
difference c*oefficient 
frequency 
channel height 

acoustic- intensity in axial direction 

yn 

total number of grid points in y direction 
total number of gprid points in x direction 
length of duct 
spinning mode number 
constant 

transverse mode number 

dimensionless Fourier coefficient of pressure p(x,y), 

2 

amplitude of pressure fluctuation or p*cj^ 

pressure defined by Eq. 51 

pressure amplitude, Eq. 10 

averaged |p| , Eq. 11 

discrete pressure at ky grid point 

radial coordinate 

dimensionless acoustic particle velocity in x direction u*p*to*H*/p^ 
dimensionless acoustic particle velocity in y direction 


weighting factor 

dimensionless axial coordinate, x*/H* 
axial coordinate 
axial grid spacing 
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y (limcnisionlcss transverse coordinate, y*/ll* 

y* transverse coordinate 

Ay transverse f^rid spacing i ' 

z* acoustic impedance 

c error in pressure 

t specific acoustic impedance 

I] dimensionless frequency, f*ll*/c* or a)*H*/27rc* 
0 parameter, Eq. 29 

\* wavelength 

p* density 


axial error variation parameter 
angular coordinate 
circular frequency 


Subscripts 


exit condition 


axial index, see Fig. 1 
trajisverse index, see Fig. 1 
cell index 

harmonic component 


entrance, x = 0 


Superscripts 

* dimensional quantity 

(“) complex conjugate 

(1) I’eal part 

(2) imaginary part 


NFMERICAI. APPROACH 


Boundary Value 1-roblem 

The calculation of the propagation of sound in ducts using numerical finite 
difference and finite clement techniques has been described in Refs. 6 to 12. 
Figure J shows a typical finite'- difference grid network used in Ref. 6 to study 






5 

the propagation of sound in a two-dimensional duct. In this case, the derivatives 
in the propagation equations are expressed in terms of the values of pressure at 
each grid point. For the boundary conditions, an entrance pressure profile, wall 
impedance, and an exit impedance are required. The collection of the various 
difference equations at each grid point form a set of simultaneous equations which 
are solved to determine the pressure at each grid point. A similar procedure 
is used in the finite element theory. In general, the finite element technique is 
more convenient to use for complex geometries (ref. 9), although it requires 
larger solution times than the finite difference techniques. 

The grid system does not have to be confined to the region inside the duct, 
as shown in Fig. 1, but It may also be extended to regions external to the duct, 
as depicted in Fig. 2(a). Grid points may be extended into the far-fleld by a 
direct placement of grid nodes or elements in the far-field, such as used in 
Ref. 12 for a flanged horn or by some type of mapping function (Ref. 7). 

However, the boundary value approach presents some severe difficulties. 
Generally, numerical solutions to the boundary value acoustic problems are 
limited to relatively low frequency sound and relatively short ducts, because at 
least 12 axial grid points or elements are required per wavelength of sound in 
order to accurately estimate the acoustic intensity (Ref. 6). A wave envelope 
analysis was developed in Refs. 8 and 10 to overcome this difficulty; however, 
the technique has only been applied to the simple cases of no flow and plug flow. 

Initial Value Problem 

In an initial value problem both the pressure and impedance either in the 
far-field or at the noise source are known. This is illustrated in Figs. 2(b) and (c). 
Using a marching technique, the pressure at the points adjacent to the boundary 
can be immediately calculated. The matrbc algebra and computer storage as- 
sociated with the boundary value problem is unnecessary for the initial value 
problem. Since no matrix of coefficients has to be stored and manipulated as 
in the boundary value problem, marching techniques are not limited to low 
frequency problems. However, associated with marching techniques is the prob- 
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lem of numerical stability; that is, preventing the growth of the round-off error’s. 

Next, the governing equation will be presented along with the formulation of 
the initial value approach. 

ANALYTIC MODEL 

Governing Equations and Boundary Conditions 
For a two dimensional duct with no flow, assuming that the pressure is a 
simple harmonic function of time and that no sources exist in the medium, 

the linearized gas-dynamic equations (Ref. 13, p. 5) of continuity, momentum, 
and energy reduce to the dimensionless Helmholtz equation 
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The sound power that leaves a duct and reaches the far-field is related to 
the instantaneous axial intensity at the duct exit. The intensity can be expressed 


in dimensionless form as 


I„ = Ueal(pu) 


where the bar denotes the complex conjugate. The total dimensionless acoustic 
power is the integral of the intensity over the transverse dimension y 


Ex= ' Ix(x,y)dy 


By definition, the sound attenuation (decrease in decibels of the acoustic power 


from X = 0 to x) can be written as 


AdB = 10 log, 

E 

o 


The pressure amplitude |p| is a quantity to indicate how the acoustic pres- 
sure varies in the duct. By definition 


Ipl = + P^^^^ 


where p^^^ and* p^”^ arc the real and imaginary parts of p respectively. A 


useful quantity in comparing the boundary value and marching solution is 
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'P'x = 


'/pp dy 


which gives a measure of the pressures across the entire transverse dimension 


of the duct. 


Finite Difference Formulation 


Instead of a continuous analytical solution for pressure, marching solutions 
determine the pressure at isolated grid [X)ints by means of finite- difference ap- 
proximations. The differential governing equation, in this case Kq. 1, will be 
expressed as an algebraic equation in terms of the pressures at the gi'id points. 
In the initial value problem to be considered here, the pressures at the grid 
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points adjacent to tlic loiown boundary are found from the known pressure and 
impedance at the boundary . 

The governing difference equation can be developed by an Integr/Ulon process 
in which the Helmholtz Eq. 1 is integrated over the various cell areas shown in 
Fig. 3 


Cell 


[ l £ + (2^)2 

a/ 


dydy = 0 


area 


( 12 ) 


The finite'- difference approximations for the various cells shown in Fig. .3 arc 
expressed in terms of the coefficients 

®m*\-l,j j-1 ^ ^m*^k, j ^ *^m^, j+i ®m*^+l,j ~ ® 

The subeript m denotes the cell number. The expressions for the coefficients 
for the various cells can be found in Table I. The derivation of cell 1 of Table 1 
is presented in Appendix A of this paper. The derivation of the coefficients for the 
other cells can be found in Ref. 10 (Appendix D) and will not be presented here. 

The explicit marching technique is based on the face that pj^_j^ j is the only 
unknown. Therefore, the pressure at (k-1, j) can be solved directly as 

'•■‘‘mPk.j-H 

The above equation can be used in marching backward from the end of the duct 
as shown in Fig. 3(a). On the other hand, if we start from x equal 0 at the 
duct inlet and move towards the duct exit, then j would become the unknown. 
In this case, the marching could bo started with cell Eqs. 7, 8, and 9 of Fig. 3(b). 

STABILITY CONSIDERATIONS 
Stability Problem 

Equation l i was applied to the problem of a hard wall duct with a plane 
wave existing at the exit (L*/H* = 1) and an exit impedance equal to - D* 

A solution for Ip! defined by Eq. 11 is displayed in Fig. A for a frequency 

?V 
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parameter of 7]= 1. The known analytical solution is given in Ref. 8 for a plane 
wave in a hard wall duct 

p = cos 2rr]x - i sin 2ht;x (15) 

or substituting into Eqs. 10 and 11 yields 

lplx=l ( 16 ) 

As seen in Fig. 4, the numerical values of |p|^ are in close agreement with 
the analytical value until the axial position, x equals 0.4, is reached. Obviously, 
an instability occurs as the solution is marched from the exit (x = 1. 0) to the en- 
trance (x= 0). Clearly, the explicit marching equation as given by £q. 14 is 
unstable in this case. Consequently, a detailed analysis of the stability of the 
Hemholtz equation is required, to determine when instabilities can occur and 
how they can be avoided. 

Von Neumann Method 

To analyze the stability of the difference approximation to the Helmholtz 
equation, we must study how the errors introduced into the computation propagate. 
We may express the exact solution to the difference equation in the form 

p(exact) = p(computed) +e (error) (17) 

The error may be due to round-off or a computational mistake. The propagation 
of errors through the computation is governed by the original difference equation. 
Assuming the initial and boundary values for the error are zero, then the equa- 
tion governing the propagation of errors is the homogeneous form of the defining 
difference equation, Eq. 13. 

For difference relations involving constant coefficients, the errors may be 
expanded in a finite Fourier series of the form 

V 1. bury. 

'k,J “) /n ® ’ <“> 

n 

This expression for the error was first used by Von Neumann (Ref. 14 and 15). 

The problem of stability is studied by noting the behavior of the coefficients 
^(n). If any ^(n) is such that 






10 


w ■ . 










l</>(n)l>| (19) 

for any n, then the corresponding error harmonic would grow beyond limit during * 


during marching. 


Stability Analysis 


The stability analysis will be performed for the two marching schemt u 
shown in Kig. 5, In these cases, the (k+1, j) values are considered the wknown, 
in contrast to the earlier examples where we were moving from the exit to the 
entrance of a duct. The explicit schemes shown in Fig, 5 allow the pressures 
to be calculated from a simple algebraic equation of the form of Eq. 13. 

To cover these two grid patterns, the stability analysis will be developed 
for the grid network si.own in Fig. 5(b), In this network, the difference approx- 
imation for the Helmholtz wave equation (see Appendix A) becomes 

Ax 2 

Pk+i,j-2Pk,j = «>k,j+i-2Pk,j + Pk,j.i) 




- pj. j 


where 




The weighting factors Wj^ and ^ indicate the importance of a particular 

grid approximation in the model. For example, if equals zero, then 

Eq. 20 would represent the Explicit x grid of Fig. 5(a). 

Since the propagation errors are also governed by Eq. 20, Eq. 18 (for a 

particular n) is substituted into Eq. 20 with €. . replacing p. .. Dividing 

inT.y. 

through the i-esulting equation by regrouping terms with the like 

k 

values of (p and noting that 




• -.yv 


e - 2 + e = -2(1 - cos nt^y, = -4 sin — 


then the solution for (p obeys the quadratic form 


a(p - 2(p + c = 0 


where 


i 1 + 0 W. /2 - d/2 ! 

i. J 


1 + e Wj ^/2 - d /2 


d = (27I7 )Ax) 

\^y/ 2 

The well known solution for tp is written as 


So that the errors will not grow, \0 | must have a value less than unity 

(see Eq. 19). This requires that the coefficients a and c obey the following 

inequalities 

-1 < - < 1 (31 

a 
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However, it is convenient to use a more restricted form of Eq. 31 
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Eqs. 32 and 33 put stability restrictions on the maximum allowable size 
of axial spacing Ax; however, Ax must be restricted to even smaller values to 
reduce truncation error. From Kef. 6, the required number of axial grid points 
was 


Thus, 


K --- Nt; Ef/li* N > 12 


K 12?) 


Also, the transverse spacing Ay was given in Ref. 6 by 


Ay = 


1 

J-1 


(34) 


(35) 


(36) 


However, no requirement for the size of transverse spacing Ay was given in 
Ref. 6 other than the number of transverse grid points be increased until con- 
vergence is achieved. 

Finally, whenever Eqs. 32 and 33 are conservatively applied, the sine 
term in Eq. 29 is assumed to take on the maximum possible value 

sin^ mr^ ^ ^ ^37) 

2 

Now, Eqs. 32 to 37 will be applied to the grid forms in Fig. 5. 

Explicit X-grid. - For the simplest explicit X* grid (Fig. 5(a)), the weighting 
factor Wj^ j is zero and Wj^ ecjuals 1. Thus 

a = ' 1 + (0/2) - (d/2)!"^ (38) 

j 

and 

e/a = 1 (39) 

Applying the first co.idition for stability, Eq. 33, yields d> 0 where 0 and d 
are defined by Eqs. 28 and 29 respectively. 
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(40) 


Finally using Kq. 3? yields the following requirement for numerical stability of 
the 5 point X-grid shown in Fig. 5(a) 

Ay>iLJ12 > V-?? = U.32 (41) 

V P/c* 

\ 

Along with Eq. 41, the second condition for stability, Eq 32, most also hold, 
but the right hand side of Eq. 32 can nrt be satisfied since c/a in Eq. 39 is 
identical to 1. In this case, the solution for </? occurs as complex pairs of 
magnitude equal to unity 

(p ” 1, if Ay > and Ax < (42) 

7) 127) 

Such a case is defined (Ref. 14. p. 337) as having a linear instability. 

At first gliuice, these results, Eqs. 41 and 42, are disappointing. Consider 
Eq. 41. With numerical techniques it desirable to have Ax and A>, related 
so that both Ax iuid Ay can both be decreased to zero to check convergence. 
However, in this special case of the I'elmholtz cxiuation, this limitation on Ay 
has an import;mt physical interpretation. Also, Eq. 42 has the exact property 
(tp - 1) required for acoustic calculations. Both Eos. 41 and 42 will now be 
considered, beginning first with Eq. 42. 

The Fourier components for propagating transverse acoustic modes should 
travel with undiminished am|)litudcs in hard wall ducts (see Eq. Bi - Appendix B 
or Ref. 16). Since could represent the x dependence of a pressure wave 
just as well as an error, tp must be held at one in the numerical calculation. 

Thus, Eq. 42 has exactly the properties desired in modeling acoustic propagation 
down a hard wall duct for propagating transverse modes. That is, the amplitude 
of the waves in the x direction is constant (tp - 1). Other nodal systems could 
be used in this problem besides the two shown in Fig. 5. However, if tp should 
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be significantly different from one, this metliod would lead to incorrect answers. 

Unfortunately, if <p equals one, errors are also not diminished In magnitude 
either. Errors at successive steps may accumulate and because of linear in- 
stability the amount of round-off error may ultimately beeome the dominant factor 
in the calculation if many grid points arc used (Ref. 15). 

Explicit U-grid . - To prevent the accumulation of errors, the 7 point or ex- 
plicit IJ-grid will be employed, as shown in Fig. 5(b). In this case, the weighting 
factor is kept at a small value (10“ ) so thet tp will be slightly less than 

one so, that errors do not grow. However, because is such a small value, 

<p can be considered equal to one for practical purposes. Consequently, the 
decrease in amplitude of the desired propagating transverse mode will be negligible. 

To first order, the stability criteria for Ax and Ay are the same as for 
the 5-point or explicit X-grid shown in Fig. 5(a). The explicit X-grid will, how- 
ever, still be used to start the marching solution. In such a case, the error as- 
sociated with a single step would be negligible. 

Now, let us eturn to Eq. 41 for a physical interpretation of its meaning. 
Assume that the argument of the sine term in Eq. 40 is small such that 

(43) 

2 2 

Substituting Eq. 43 into Eq. 40 yields 



r 

tc 




mu’ 


n < 2fj 


(44) 

KF* ■ 



Eq. 44 represents the cutoff condition for a mode to propagate unattenuated down 
a hard wall rectangular duct (see Eq. B2, Appendix B). Therefore, Eqs. 41 and 42 
can be Interpreted to mean that only propagating modes can be handled in the 
calculation. 

Using Eq. 36 with Eq. 41, the Increment Ay can be expressed in terms of 
the number of grid points J in the y direction as 

J < TT() + 1 (45) 

For a given harmonic to exist in the solution, there must be sufficient grid points 
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to resolve a Fourier component. To define a wave, at least 2 points arc re- 


quired. Therefore 


-V.J 

n < — 


or substituting Eq. 46 into Eq. 45 yields 


mj + 1 


ft 


P*' 

m 


which is roughly equivalent to Eq. 44 for small values of >j. Thus, the stability 
criteria indicate that marching c-ould be imstable if sufficient grid points are 
available to resolve the higher order transverse modes. 

As will be shown later, the restriction of handling only propagating modes 
will limit the useful application of this technique to high frequency (hi^ rj) 
calculations and to those problems where the cut-off modes do not play an 
important part in matching pressure or other boundary conditions. 

DIFFERENCE ALGORITHM 

The marching equation for calculating pressure down the duct for the central 
cell caji now be written from Eq. 20 as 


.-y £'■ 




Pk+i,j <2 - (27n;Ax)2) p,^ . - P|^_j j (-PkJ+1 ' Pk,j-P 




mf 

sm I if' 




For the upper wall cell, soft wall boundary conditions require that 


PkH.i - Pk, j - Pk-ij ^ '"k ^ ^ »k,i * 


o h<,j ^ Pk,j “ hc,j-l 


+ 2 W., 


Ax i27n)Ay 


j Pk-ij '■Pk-i.j- Pk-i.j-i 

*'k-l 





I - ^ 
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The fuation at the lower wall requires that j - 1 in Eq. 49 becomes j + 1. For 
—6 

Wj^_^ equal to 10” , the equations represent the explicit H-grid system. For 
Wj^ l equal zero and Wj^ equal unity, the equations represent the explicit X-grid. 

DUCT MARCHING EXAMPLES 


To illustrate the capabilities and limitations of the marching technique, 
sample problems will be presented for the simplest case of a two-dimensional 
duct as shown in Fig. 6, instead of the more complicated geomtry shown in Fig. 2. 
The advantage of this simplification is that a direct comparison of numerical and 
analytical techniques can be made. A boundary value solution for Fig. 6(a) can 
be used to set up the exact exit conditions for the problems in Figs. 6(b) and (c). 

As a result, the numerical accuracy of the marching technique can be conveniently 
checked for the case of soft wall ducts. 

Hard Wall Ducts 

In the example shown in Fig. 7, the stability condition for a hard wail duct 
is evaluated for decreasing values of Ay (increasing J) for various modal inputs 
at a dimensionless frequency t) equals 3. As seen in Fig. 7, for the first mode 
(n = 0), the calculations are convergent to the analytical solution for J considerably 
larger than the value of 10 predicted by Eq. 45. It is likely that the round off 
error is quite small for a plane wave and the calculation is terminated before the 
error can grow to a significant value. On the other hand, the higher order trans- 
verse modes (n > 0) follow the stability criteria quite closely. 

As an additional check shown in Fig. 8, the number of grid points in the 
x direction was increased to test for stability. As theory predicted, the results 
remain converged. 


Soft Wall Duct 


Next, the marching technique is applied to a soft wall duct. From Ref. 10, 

the boundary value finite difference solution for a plane wave input was used to 

establish the exit conditions. The results arc shown in Fig. 9. As can be seen 

in Fig. 9, the explicit H-grid marching technique and the boundary value technique 
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are in exact aKreemcni.. A(>ain there was no error growth when the number of 
axial grid points was increased. 

l^asc Uclationship 

A marching theory was devised so that the acoustic pressure distribution 
at a source might be deduced from far-field measurements. In actual practice, 
Refs. 17 and 18 for example, sound pressure level measurements are recorded 
at a set of microphone stations in the far-field without regard to the phase 
relationship between the various microphone measurements. In the next example, 

i 

the sensitivity of the duct attenuation and source profiles to a change in phase at 

' I 

the duct exit is considered. 

. For the same problem considered in Fig. 9, again the finite difference solu- 
tion from Ref. 10 was used to establish the exit condition. The solid line shown 
in Figs 10 and 11 indicates results based on the correct phase relationship as 
computed from the boundary value analysis. On the other hand, the dotted line 
in Figs. 10 and 11 presents a calculation with an exit pressure of the same 
absolute pressure magnitude (8PL) but of different phase. In this case, the 
imaginary part of the acoustic pressure p^^^ was assumed to be zero. As seen 
in Fig. 10, the improper phase relationship leads to much higher attenuations. 

Fig. 11 shows a pronounced effect of the exit phase on the calculated source 
pressure profile. Obviously a considerably different pressure profile from 
that of a plane wave has been generated at the liner entrance. Clearly, if one 
wishes to establish the structure of the noise source from far-field measurements, 
the far-field measurement must include a reference microphone to establish 
the phase, 

FAR-FIELD APPLICATION 

As mentioned in the introduction, the impetus to develop the marching tech 
technique came from a desire to establish the noise source from far-field data. 
The purpose of this section is to discuss the possible application of this techni- 
que in obtaining information on the spinning mode content of a turbofan engine. 

In analyzing far-field data from a typical turbofan engine, cylindrical and 




spherical type coordinate systems should bo used to resolve the fan's radial, 
circumferential (inside duct) and axial pressure. In such a case, a three 
dimensional differential wave equation of the form (inside duct) 

L(p(r,^(»,x)} = 0 (50) 

is required. The differential operator L should, of course, be sufficiently 
general to account for variable mean flow Mach numbers and variable inlet 
geometry. The operator could be cylindrical coordinates inside the duct and 
spherical coordinates outside the duct. The coordinate system would be 
coupled at the duct exit using some type of matching (see Ref. 19). 

Inside the duct, the circumferential variable ^ could be separated by 
assuming a pressure distribution of the form 

im 

p(r, i^,x) = Pj^(r,x)e ® (51) 

where m is the "spinning" mode number. Therefore, the equation describing 
s 

the pressure field reduces to a simplier two-dimensional form, 

L^p^(r,x)}=0 (52) 

Now, let us consider the possible application of Eq. 52 to a numerical solution 
for the noise source based on far- field data. 

To estimate the modal content of a fan source using the marching technique, 
a matrix of grid points would be placed in the acoustic field leading from the 
microphones to the fan source. Assume that the microphones are spaced about 
a fan inlet at some fixed height above the ground under free field conditions. 
Furthermore, assume the microphones are set to measure both magnitude and 
phase of the pressure. As was shown in the last example of the previous section, 
phase information is necessary to account for cancellation of the various modes 
in the duct and the far-field. 

Unfortunately, a unique solution exists at the fan source for any assumed 

mode number m^; that is, a lobe pattern in the far-field would transfer (by 
s 

marching) to a radial shape in the duct at the source. For numerical calculations. 
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the ranse of m would be limited to low values of m (Ref. 3) so that some 

a S 

propagating radial modes exist in the duct. For large values of m , most if 

s 

not all of the radial modes no longer propagate as a wave (Ref. 16, p. 497). 
Although analytical solutions might allow predictions of the noise source for 
large m^, such a solution would be no practical interest since the pressure 
would propagate by diffusion (like heat) and errors in measurement would be 
highly magnified (Ref. 20, pp. 229-231). The Von Neumann analysis will 
again limit the grid spacing to allow only propagating modes in the duct. 

Although only propagating modes exist in a duet near its exit, higher order 
far-field modes can be generated at the exit by a change of geometry (Ref. 21, 
p. 122 or Ref. 22), in this ease the inlet lip of a turbofan engine. Therefore, 
since all modes propagate in the far-field (Ref. 13, p. 212), in applying the 
marching technique to the far-field, the grid structure must account for a full 
raiige of acoustic modes. 'I'his will be the case since the number of grid points 
can increase because of the inlarged area available in the far-field. That is, 
using a rectangular system for example. Ay is constant both inside and outside 
the duct, Eq. 42; consequently, more grid points can be used in the far-field. 

Because the spinning mode number m„ is not known, and the source will 
be multimodal, measurements from a far-field horizontal ari'ay of microphones 
alone can not uniquely determine the circumferential modal content of a fan. 
bi effect, a horizontal traverse in the far-field is analogous to a radial traverse 
internal to the duct, neither of which give any information on the circumferential 
modal content of the source. This agrees with the theory of Ref. 5. 

If the marching method is to be used to deduce circumferential modal patterns 
inside the duct without using any apriori assumptions about the modal content, 
additional pressure measurements must be taken in the far-field in the colatitude 
direction (Ref. 23, p. 60). A three dimensional solution of the complete set of 
far-field data would lead to a source distribution p(r, 0) at the fan face. This 
would be decomposed (Fourier scries) into the various circumferential modes. 

If the phasing among the modes is a random variable such us for broadband 




noise or if a large number of modes are involved even for a discrete tone, the 
direct solution for duct modes from far-field measurements might not be practical. 
For these cases an approximate approach such as used in Ref. 3 or 4 mig^it be 
used. 

SUPPRESSOR ATTENUATION WITH PLANE WAVE INPUT 

The noise attenuation of a suppressor is now evaluated using the marching 
techniques. Consider the attenuation produced in the beginning length L* of an 
infinitely long duct which has a uniform impedance along its entire length, as 
shown in Fig. 12(a). To numerically predict the attenuation of the entrance 
length L* in the infinite duct (Ref . .6 to 10), a finite duct of length L* (Fig.. 12(b)) 
is used with an assumed exit impedance hopefully equal to the impedance which 
exists at L* in the infinite duct. 

An infinitely long duct with uniform impedance will not have any reflections 
at any position in the duct. ConscKiuently, the wave propagation in the entrance 
region of the infinite duct can be represented with a numerical solution to a finite 
duct of length L*, as in Fig. 12(b), by choosing the exit impedance at L* so 
that the reflections are minimized. For the case of zero mean flow, and a plane 
wave input, a p*c* exit impedance was found to give good agreement between 
the analytical and numerical results for 77 between 1 and 5 and L*/H* between 
0.5 and 6 (Ref. 8 , Fig. 9). 

Herein, the axial impedance of the medium is assumed to be constant down 
the duct; consequently, the impedance of the medium at the noise source (x = 0 ) 
is also assumed to be PqC* . This entrance impedance assumption precludes any 
possibility of reflections at the duct entrance and would not be valid in cases 
where strong reflections are known to exist. Using this assumption and the 
marching technique, the noise attenuation at the optimum point (point of maxi- 
mum attenuation in the impedance plane) is now calculated for a two-dimensional 
duct with L*/H* values ranging between 0. 5 and 5 and an input plane wave with 
dimensionless frequencies rj of 1 , 2 , and 5. The values of the optimum liner 
impedance used in the numerical analysis were determined from the analytical 
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tofhniquc's prescnlod in !{ef. 2-1. 

Tne numerically caltulaled attenuations are compared to the corresponding 
analytical results in Fig. 13. The numeric.nl results are given by symbols while 
the analytical results are represented by the solid lines. For dimensionless 
frequencies t] equal to 1 or 2, the numerical and analytical results are generally 
in poor agreement. Also, for i] equal to 1, the finalysis breaks down (calculated 
numbers go to infinity) for L*/l\* equal or greater than one, while for rj* equal 
to 2 the break down begins at L*/ll* c*qual to .3. Most likely, the truncation 
error is large because there are few points in the transverse y direction. In 
the case of 7) equals 1, Eq. “15 allows the use of only 4 points. 

On the other hand, the analytical and numerical results are in very good 
agreement for 7) equal to 5. It is well known (Ref. 25, p. .391), that the trunca- 
tion error is dominated by the largest spacing increment, in this case Ay. For 
high frequencies, sufficient transverse grid points are now available to reduce 
the truncation error in the y-direction and thereby to accurately predict the 
transverse pressure profiles. Therefore, stability restrictions on Ay limit 
the application of the marching technique to high frequency sound, r\ equal to 5 
or greater. 3'his limitation applies equally well to the work of the previous 
section (marching from the duct exit to the entrance). 

In the simple but practical example just considered, reflections in the duct 
do not exist. In such a case, the aerodynamic forces producing the noise source 
could be calcul.nted without reference to the acoustic field. However, a wide 
class of problems exists in which the acoustic field reacts back in the source 
(Ref. 2G). such as when partial blockage (Sc.attering) occurs in a duct (Ref. 16, 
p. -40U) or when a duct mode produces a large back reaction force on the source 
(Ref. 16, p. 501). In these more complicated situations, an iterative marching 
approach might account for reflections in the system. 
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CONCLUDING REMARKS 

1. An initial value numerical method is developed to predict the pressure 
distribution of a ducted source from the downstream pressure measurements. 
The technique is useful for relatively high frequencies ( 7 / s 5). For the special 
case where reflection in a duct are small; tlie technique can also be used to esti- 
mate the attenuation of a straight soft wall duct. 

2. Compared to conventional finite difference or finite element boundary 
value approaches, the numerical marching technique is orders of magnitude 
shorter in computational time and computer storage requirements. Therefore, 
this technique could be used as a companion to the boundary value techniques 
particularly at hi^ frequencies. 

3. The Von Neumann method accurately predicts the onset of instability in 
the numerical calculation. The stability considerations limit the approach to 
high frequency sound. 

4. Because the marching technique is explicit in nature, the technique 

can readily be applied to three dimensional problem. Thus, the techni<"ie could 
be useful in tracing.far- field pressure patterns back to the noise source. Ex- 
ample problems also indicate that the phase relationships between microphone 
measurements in the far- field must be recorded if the source distribution is 
to be calculated. 

5. The marching technique, in effect, reduces the partial difference equation 
to an ordinary differential equation because the transverse gradient terms are 
approximated by previously known values. Thus, in principle, a nonlinear prob- 
lem could be conveniently handled with Runga Kutta or other standard techniques. 
In addition, the techniques may be applicable to other complications in the flow 
field, such as a variable speed of sound (Ref. 27). Of course, the procedure 
would have to be first applied before one could be certain of its success. In 
these more complicated cases, the stability analysis could be performed by 
computer trial and error to determine critical values of Ay and Ax. 


2 :? 





6. The marching technique can easily be applied directly to the Individual 
continuity, momentum, and energy equations. Thus, it is unnecessary to com- 
bine these equations to obtain a wave equation. The combination was performed 
herein to simplify the stability analysis. Hopefully, with the addition of flow 
gradients and some nonlinear effects, the stability requirements will not be as 
restrictive. 

7. Finally, the marching technique can not be used directly to estimate the 
attenuation of noise suppressors in situations where large reflections occur in a 
duct such as would occur with a converging duct. In this case, an assumed 
entrance impedance of p*c* could be invalid. However, an iteration process 
might be developed which could lead to a valid solution. This limitation does 

not apply to situations where both the pressure and entrance impedance are known 
from measurements. 
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APPENDIX A 

FINITE DIFFERENCE FORMULATION 
The governing difference equation can be developed by an integration process 
in which the Helmholtz Eq. 1 is integrated over the cell area, Eq. 12. Eq. 12 can 
be rewritten as 


dxdy + 

2 


dx 


Cell 

area 


^ dxdy 


9y 


(23rrj)^ pdxdy = 0 


(Al) 


Cell 

area 


Cell 

area 
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The various pressures and their derivatives in Eq. Al are evaluated at the edges 
of the cells, designated by + and - signs in the following equations. These terms 
will be different for each cell because the boundary conditions, as well as size 
of the cell vary. 

For cell number 1 in Fig. 3, the standard approximations for the derivatives 
are used 


I ^ riydv = 




' / ^dxdy 

I 9 

. I dx^ 


dx 


' '+ 
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' pdxdy = (27r77)“ Pk j ' " (27rrj)^ AxAy pj^ ^ (A4) 


where 
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' _ Pk+l. j ~ Pk. i . 
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(A5) 
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^.i+1 ~ i . ^ - PkJ-l 

dy Ay 9y Ay 


Substituting Eqs. A2 to A5 into Eq. A1 gives 


^ \ 1 --Pi. • ,+2'i+;J^\<2222M! 
Ax lAx 2 




The finito- difference approximation for cell no. 1, Eq. A7, as well as the other 
cells in Fig. 3 can be written in the form of Eq. 13 in the body of this report. 

The derivation of the difference equations at the upper boundary are just as 
easy provided one recognizes that 

<>y f,„ 


as given by Eq. 5. Similar expression can be developed for the exit condition. 

For the special case of the explicit Il-grid shown in Fig. 5(b), the transverse 
grid derivatives at the upper wall would be approximated by 

ilfil = ^Ic^Pkj ~ PkJ-lj ^k-l^Pk-lJ~ Pk-lJ-1^ . . j 


^- -Wki2T7Pk i W^^_^i2 ,;p,^_^ 

<^y >-1 y \ 


(AlO) 


The weighting factors Wj^ and Wj^ are defined by Eq. 21 in the body of this 


paper. 
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APPENDIX B 

PIlYblC'AL SICNIEICANCE OP VON NEUMANN METHOD 
The Von Neumann method Ih quite compatible with acoustic problems since 
the Fourier representation of the error at each grid point can be also applied to 
an actual transverse harmonic pressure input at the boundary. Therefore, in 
addition to evaluating the stability of the various grid systems, let us consider 
what effects the governing equations and boundary conditions will have on a 
transverse harmonic input for the simple case of a hard wall duct. 

The solution to the Helmholtz equation, Eq. 1, for a two-dimensional in- 
finitely long duct with rigid walls is given for example in Ref. 16 (pp. 496 and 
504). The nondimensional acoustic pressure can be expressed (using notation 
different from Ref. 16) as 
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P = 


cos nTry e 


-iff '/(2 tj) 


2 

n y 


(Bi) 


n=l 


where the rigid walls are located at y equals 0 and 1 . The symbol n represents 
the transverse mode number. The characteristic function cosnyy satisfies the 
condition that the transverse velocity at both walls will be identical to zero for 
a hard wall two-dimensional duct. 

Eq. Bl indicates that for the mode number n, if 


n < 2j) 


(B2) 


the mode (cosiiffy) will propagate unattenuated down the duct with a constant am- 
plitude. On the other hiuid, if n is greater than 2rj the mode will exponentially 
decay. It is said to be cutoff. For the special case of n equals 2 t) the nr.ode 
becomes a transverse standiu ; wave in the duct; that is, no variation exists in 
the X direction. 

In a hard wall duct, if Eq. B2 is satisfied, chc numerical marching technique 
should predict that the amplitude of a transverse harmonic wave will remain 
constant along the duct. This is precisely what occurs when Eq. 42 is satisfied; 
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that is, 0 m':st b*; equal to one in the Von Neumann analysis. If 0 were less thaii 
one, the transverse harmonic mode would decay with distance, which would violate 
the true physical situation as indicated by Eq. B1 for propating transverse pressure 
modes. Thus, the numerical marching technique has the proper characteristics to 
match the physical situation only if <p equals 1 in the Von Neumann analysis. 

In the sample problems investigated in the body of this paper, plane wave 
and single transverse mode inputs were considered. I-'or a plane wave input. 

Ref. 6 had found that an exit impedance of is a good choice for both hard 

and soft wall ducts. However, for a single harmonic mode, the exit impedance 
associated with Eq. B1 is used. Substituting Eq. B1 into Eq, 2 yield for the axial 
acoustic velocity of the nth mode. 

^ (B.3) 


E’inally, substituting Eqs. B1 ;md B3 into Eq. 6 yields for the exit impedance 

(B4) 

For a single harmonic mode n, the use of Eq. B-i as the exit impedance will 
allow numerical calculations in finite ducts to be compared to the analytical re- 
sults in an infinite duct. 
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rSPECiFIC ACOUSTIC WALL tMTOMIVeC. C 
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Figure 1. - Grid-point repi mentation of tv^-dlmenslonal flow duct. 
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Figure Z ■ Types of boundary and Initial wlua problems 
Ibr acoust'c pressure fields. 
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AVERAGE DUCT ACOUSTIC PRESSURE. |p|. 


rCiil NO. 2 .'•WAU BOUNDARY 



(a) BACKWARD MARCHING. 
WAa BOUNDARY 



Figure 3. * Coordinate and explicit grid point representation of 
of acoustic field. 
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Figure 4. • Averaged integrated pressure distribu 
tion |Plx in a hard wall two dimensional duct 
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Figure 5. - Explicit difference relationships for the Helmholtz 
equation. 


» plane 

t 

V« : : : : : : 1 

lal BOUNDARY VALUE 
PROBIEM. 

• ••• 

- _ - * «• * 

• ; • ; X 1 1 1" ’"o"® 

<)| FAR-HELD INITIAL 
VALUE PROBIEM. 

v: I : 

m . > . 


(Cl SOURCE INITIAL 
VALUE PROBIEM. 


Hgure 6. - Propagation problems tor a ducted 
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Figure 7. - Stability evaluation in marching from 
exit to entrance for hard wall duct It? ■ 3, 

1. K - 36). 
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Figures. ’ Stability evaluation in 
marching from exit to entrance 
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Figure 11. - Pressure profiles at duct entrance 
as a function of the phase at exit. 
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(a) ANALYTICAL MODEL 
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(b) REPRESENTATION OF NUMERICAL EQUIVALENT. 

Figure 12. - Numerical representation of the attenuation problem 
for a section L* of an Infinite duct. 
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Figure 13. * Effect of axial length and 
frequency on attenuation at opti- 
mum impedance in two-dimensional 
duct for a plane wave input. 
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